Phase transition in a 2-dimensional Heisenberg model 
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We investigate the two-dimensional classical Heisenberg model with a nonlinear nearest-neighbor 
interaction V{s,s') = 2K[{1 + s- s')/2]^. The analogous nonlinear interaction for the XY model 
was introduced by Domany, Schick, and Swendsen, who find that for large p the Kosterlitz-Thouless 
transition is preempted by a first-order transition. Here we show that, whereas the standard (p = 1) 
Heisenberg model has no phase transition, for large enough p a first-order transition appears. Both 
phases have only short range order, but with a correlation length that jumps at the transition. 
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The two-dimensional Heisenberg and XY model are such close relatives that it has taken a long history of efforts 
before their properties could be told apart. Both are special cases, for n — 3 and 2, respectively, of the 0(n) symmetric 
Hamiltonian 

H = -KY,^^- (1) 

<i,j> 

Here Si is an n-component spin of unit length at lattice site i, the sum is on all pairs of nearest-neighbor sites of a 
two-dimensional lattice, and K = J /k-oT. For all n > 1 the system has c? = 2 as its lower critical dimension. 

Bloch's 1930 spin wave argument Q], put on a firm mathematical basis only much later by Mermin and Wagner 
implies that neither the XY nor the Heisenberg model can have a spontaneously magnetized low-T phase. The 
early investigations dealt exclusively with the Heisenberg model. In 1958 Rushbrooke and Wood, after studying high- 
T series [Q, first remarked that in spite of Bloch's argument the possibility of a phase transition in the Heisenberg 
model should be taken seriously. This was reemphasized in 1966 by Stanley and Kaplan j5|, who envisage, for the 
Heisenberg model, a low-T phase with an infinite susceptibility 

In the late 1960's the high-T series of the Heisenberg the XY model were compared [^,0. Qualitative similarity 
was found, but no general agreement was ever reached about the significance of certain quantitative differences. A 
phase transition in either model continued to be considered by many as only a remote possibility, until Kosterlitz 
and Thouless (KT) |^ demonstrated that there is a phase transition in the XY model and clarified its topological 
character. 

Since the KT arguments were specific for n — 2, the two-dimensional Heisenberg model (and, indeed, the Hamil- 
tonian (|^) for all n > 2) has from then on been believed to be without a transition. Further support for this view 
came from the analytical low-T renormalization group approach developed by Polyakov Brezin and Zinn- Justin 



10[ , and Nelson and Pelcovitz 1 11 1 , and from Monte Carlo renormalization due to Shenker and Tobochnik |_2[ . The 
absence of a rigorous proof has however left room for arguments ( [|l3| and references therein) that the Heisenberg 
model (Eq. (|l]) with n = 3) may after all have a phase transition; this is not, however, our point of view. 

Here we consider the 0(3) symmetric Hamiltonian 

^ = - E ^(^"^ • (2) 

<i,'J> 

where V is an arbitrary nonlinear function. For reasonable choices of V (in a sense not a priori clear) one expects 
that (H) is in the same universality class as the standard "linear" 0(n) model (|l|). Expression (||) is interesting for at 
least two reasons. 

First, the freedom to choose F is a key ingredient in theoretical analyses by Villain |l4j of the 0(2) model and 
by Domany et al. [|l^ and Nienhuis of the 0{n) loop model. For n > 2 the latter model does undergo a phase 
transition [|l7| which corresponds to a hard-hexagon-like ordering of the loops. But in spin language the transition 



1 



appears to occur in an unphysical parameter region with negative Boltzmann weights. It does not provide evidence 
for a phase transition in 0(n) spin models with n > 2. 

The second reason of interest in (0) comes from the relevance of the KT theory for the melting of thin adsorbed layers. 



The difficulty encountered in observing the predicted |18| hexatic phase, whether experimentally or in simulations, 



was suspected by some to be due to the KT transition being preempted by a first order transition as a consequence 
of various nonlinearities not incorporated in the theory. Domany et al. (DSS) therefore investigated an 0(2) 
symmetric XY model with a specific nonlinearity controlled by a parameter p, viz. 

V{s,-s,)^2K[il + s,-Sj)/2]P (3) 

(our p is their p^). Indeed DSS found by Monte Carlo simulations that for strong enough nonlinearity {p ^ 50) the 
KT transition is replaced with a first-order one from the massless low-T phase to a high-T phase with exponentially 
decaying correlations. While this suggests that melting via a hexatic phase may similarly be preempted by a first-order 
transition, the DSS result has been subject to controversy |^ . 

Here we confront again the XY and Heisenberg model. We have Monte Carlo simulated the latter with the nonlinear 
interaction (||) on square L x L periodic lattices. Randomly chosen orientations are accepted with Metropolis-type 
probabilities. Slow relaxation at low T limits the largest system size to about L = 200. 

No signs of a phase transition were seen for p ?s 1, but for p = 20 there is a clear jump in the energy as a function 
of K. Fig. |l] shows the resulting hysteresis for a system of size L = 48. For the XY model a similar narrow hysteresis 
loop was observed by DSS, but today's computers yield a clearer picture in the Heisenberg case. 

Similar Monte Carlo runs for p < 20 show a weaker first-order character, but do not clearly show where the first- 
order line ends. In order to answer this question, we have determined the specific heat for a grid of points in the K-p 
plane. We thus found the specific-heat maxima as a function of K. Fig. || displays these maxima Cinax(j^;-^) versus 
L. In the absence of a phase transition Cniax(p, L) ~ est when L increases; this behavior is seen for small p. In its 
presence we expect, at large L, 

C^,^{p,L)^coL^y-^ (4) 

with y = 2 (y < 2) in the case of a first-order (continuous) transition. The data for p = 20 in Fig.|^ are consistent 
with y = 2. The finite-size divergence weakens for p < 20, and the p = 16 data indicate a continuous transition with 
y = 1.84 ± 0.05. The downward trend at even smaller p is consistent with Cinax(?3, L) ~ est at large L. This suggests 
that the first-order line in the p-K diagram ends in a critical point near p = 16. 

Simulations for p > 20 show an enhanced first-order character. Transition points were found by several runs, 
starting with half the system fully aligned, and the other half chosen randomly. The results, which hardly depend on 
L for L > 32, are shown in Fig. |^ versus p. 

The transition points can also be estimated from the high- and low-T expansions of the free energy. Neglecting 
loop diagrams in the high-T expansion the lattice effectively reduces to the Bethe lattice (BL). Its partition function 
'per bond' is 



/oo 
ds exp 2K[{1 + s-i)/2]P = 47r^ — 
I — n V 



(2i^)* 



(5) 



where the prefactor accounts for the phase space volume of a spin and the sum for the spin-spin interaction. For N 
spins and zN/2 bonds we thus have 

which yields the high-T approximation Tht of the free energy of a square lattice (z = 4) of A'' = sites as 

= - log(47r) - 2 log -^^^^ (7) 

At low T, the spin-wave approximation (SWA) of H is 

Hswa({sJ) = -^NK + Hoiisn) + nG{{sU) (8) 
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where is the Gaussian Hamiltonian 

HG{{s^}) = lpK^{s,-SJf (9) 

<v> 

By standard methods one obtains from it the low-T approximation Fjjy to the free energy, 

^ -iK - log(47r) + log{8pK) 

L-l 
m,n— 

where the prime indicates that (m, n) = (0, 0) is exchided from the sum. For large N — L"^ the sum on m and n tends 
towards —2 log 2 + 4G'/7r = —0.2200507 ■ • ■ where G is Catalan's constant. 

The intersection of the two free-energy branches was found numerically for several p. The resulting approximation 
of the first-order line, shown in Fig.|^, is in a good qualitative agreement with the Monte Carlo results. 

Next, we check the consistency of our magnetization data for the low-T phase with the Mermin- Wagner theorem 
1^,1). Fig.^ shows that the mean square magnetization m? = L^^^^^^jisi ■ Sj) decays slowly with L. In contrast, 
the energy rapidly tends to a constant with increasing L. 

In order to compare this magnetization behavior to theory, we recall that in the standard {p = I) Heisenberg model 
the correlation length ^ is well fitted |l|] at low T by ^{K) « Cexp(27ri^)/(l + 27r/v) with C « 0.01. For K r < ^ one 
expects the SWA result g{r) = {si ■ Si+r) ~ to hold, where ry — XjixK. Consequently w? ^ L^^ dr rg{r) ^ L^'' 
for 1 ^ L < For ^ < i the integral on r converges at the upper limit and one has L^^. 

Now take p ^ 1 in the model under study. Then the angle 6 between two neighboring spins is in a narrow two- 
dimensional harmonic potential well as long as <C np^^^^. For 7rp^^/^ < 6 the Boltzmann weight is decreased by a 
factor exp(— 2i4') and almost independent of 9. When K ^ 1, most angles are small, and g{r) will behave according 
to the SWA, but with an exponent 77 = l/irpK; and the correlation length £,{pK) estimated as above will exceed any 
system size L attainable in simulations (disregarding a renormalization effect of ^ due to the nonlinearity of V) . 

Next let iiT ~ 1 while still pK ^ 1. Then the fraction of nearest neighbor spins with large relative angles will no 
longer be exponentially small in K. This will cause a downward renormalization of the effective coupling of the SWA, 
if this concepts remains at all applicable, and of ^, but it is not a priori clear if ^ will still exceed the system size. To 
answer this question we consider Fig.^. For p — 20 and K = 1.4 the unrenormalized SWA gives ry = l/irpK = 0.012. 
Fig. ^ confirms the power law decay of m^, but yields a renormalized exponent rjcs ~ 0.030, estimated from the range 
32 < L < 192. This corresponds to an effective SWA coupling K^ff ~ 10.6. We note that £^{Kce) is still very much 
larger than our L values, which indicates the self-consistency of the renormalized SWA. Hence we conclude that the 
low-T phase has a correlation length ^ much larger than the system sizes L considered here, and has a pair correlation 
that, at these distances, decays as a power law. 

Our finite sizes L restrict the spin waves to small deviations, so that is considerable. One may ask how stable 
the first-order transition is under large deviations occurring in large systems. We have imposed large-amplitude 
waves using antiperiodic boundaries in both directions. This reduces considerably in finite systems at low-T, and 
renders the low-T phase less stable. Monte Carlo data a.t p — 20, i = 48 show that the energy jump and hysteresis 
are strongly suppressed. The deformation energy per bond is oc L^^. Fig. |5| shows that for L = 192 indeed the 
first-order character is partly restored to the situation of Fig. ^ This indicates that the first-order transition persists 
even when spin waves suppress the magnetization at large L. 

It is not clear how to define an order parameter reflecting a symmetry of the model. The phases separated by the 
first-order line have different degrees of short-range order, as is the case in a gas-liquid system. Thus we expect the 
first order line to end in an Ising-like critical point. Indeed, our result y — 1.84 ± 0.05 agrees well with the Ising 
magnetic exponent = 15/8. We note that the energy fluctuations of this model correspond with the Ising magnetic 
scaling field, because it is the energy that has a discontinuity at the first-order line. 

In conclusion, we have investigated a Heisenberg model with interactions that depend nonlinearly on the spin 
products. For strong enough nonlinearity there appears a phase transition. This transition is unrelated to earlier 
claims |lj], which applied to the linear case. But it does seem related to the DSS transition in the XY model in the 
following way. Adding a term jj^ki^k)'^ ^1- @ leads to crossover to the 0(2) model as 7 varies from to 00. In 
the jp plane we expect a line p = Pc{l) (with Pc(0) ~ 16) above which the transition is first order and below which 
it is of the KT type when 7 > 0. 
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FIG. 1. Energy per spin versus coupling for a periodic system of size L = 48 and p — 20. The energy discontinuity and the 
hysteresis indicate a first-order phase transition. Each data point results from lO*" Monte Carlo steps per site. The statistical 
errors are smaller than the symbol size. Jumps in the energy (see arrows) occurred while taking the data points on the vertical 
lines. 
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FIG. 2. Specific-heat versus system size L on logarithmic scales. The curves serve only to guide the eye. The 

data points apply to p = 6 (+), p = 7 (x), p = 8 (>f), p = 10 (□), p = 12 (■), p = 14 (o), p = 16 (•), p = 18 (A), and 
p = 20 (a). These data suggest that the critical point at the end of the first-order line lies near p = 16. Each data point was 

determined from several Monte Carlo runs which, because of slow relaxation, had to be long (up to about 10* updates per site 
each). This is where the bulk of the computational effort went. The errors do not exceed the symbol size. 
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FIG. 3. Phase diagram of the present 0(3) model in the p vs. K plane. The full curve is obtained by equating the high-T 
and the low-T expansions of the free energy. The data points represent our numerical results. 
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FIG. 4. Magnetization squared and energy(x) per spin versus system size L in the low-T phase at if = 1.4, p = 20. 
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FIG. 5. Energy per spin versus coupling for a system of size L = 192 and p = 20 with antiperiodic boundary conditions. The 
first-order character is still appaxent under these conditions. Each data point represents a simulation of 2 x 10® Monte Carlo 
steps per site. The statistical errors are comparable to the symbol size. 
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